Self-similarity of the Third Type in the Strong Explosion Problem 
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Propagation of a blast wave due to strong explosion in the center of a power- law-density (p oc r~ a ) 
spherically symmetric atmosphere is studied. For adiabatic index of 5/3, the solution was known to 
be self-similar, (of type I) for a < 3, self-similar (of type II) for a > 3.26, and unknown in between. 
We find a self-similar solution for 3 < a < 3.26, and give a (tentative) numerical proof that this 
solution is indeed an asymptotic of the strong explosion. This self-similar solution is neither of type 
I (dimensional analysis does not work), nor of type II (the index of the solution is known without 
solving an eigenvalue problem). 



Introduction — We consider propagation of a blast 
wave due to strong explosion in the center of a power- 
law-density, p oc r~ Q , spherically symmetric atmosphere. 
This problem is of interest for two reasons. 

An exploding star sends a blast wave into the stellar, 
and later into the circumstellar medium. The (circum) 
stellar medium has a varying but spherically symmetrical 
density, and power-law density profile would often be a 
reasonable approximation. Blast waves in a circumstellar 
medium are interesting astronomical sources (eg. [1] and 
later papers by the same author). 

The strong explosion problem is also of methodological 
interest, being a canonical example of self-similarity in 
gas dynamics. Self-similar solutions, with shock radius 
R oc t@ , of type I and of type II were known in this 
problem. 

Type I: For a < 3, the blast wave is decelerating, 
(3 < 1 [2]. The self-similar solution is of type I, meaning 
that the index (3 is determined by the energy conservation 
law - from R 3 p(R)(dR/dt) 2 oc const, one gets (3 = 2/(5 — 
a). This is a straightforward generalization of the Sedov- 
Taylor-von Neumann solution describing the explosion in 
a uniform atmosphere [3-6]. 

Type II: For a > 3.26 *, the blast wave is accelerating, 
j3 > 1 [7] . The self-similar solution is of type II, meaning 
that the index j3 is determined by an eigenvalue problem 
- it must be chosen so as to remove a singularity of the 
system of ordinary differential equations describing the 
self-similar solution. Dimensional analysis is insufficient 
in this case, and (3 should be determined numerically. 

Type III: For 3 < a < 3.26, the blast has a constant 
velocity, j3 = 1 [this paper]. The self-similar solution 
is of a new type, type III, meaning that the index (3 is 
determined by the requirement that the non-self-similar 



part of the blast wave serve as a piston for the self-similar 
part of the blast wave. Dimensional analysis does not 
work in this case, but one can calculate (3 without solving 
an eigenvalue problem, it is just (3—1. 

We consider only a spherically symmetrical model. In 
reality, deviations from spherical symmetry should exist 
[8], but are expected to be very small generically (we 
checked this numerically), should be especially small for 
a w 3, and may be entirely absent for a > 3. Gravity is 
neglected. We first calculate the new self-similar solution, 
and then describe numerical simulations which confirm 
that our self-similar solution is indeed an asymptotic of 
the strong explosion. 

Formulation of the problem 

At t = and r — » oo, the atmosphere is spheri- 
cally symmetrical, cold, motionless, and self-similar with 
p — > r~ a . Also at t = but at small r there is a non- 
zero energy density in arbitrary form. The problem is to 
calculate the resulting flow at t — > oo. 

Self-similar solution 

One assumes that at large t, the flow is a self-similar 
blast wave (meaning a shock wave followed by a flow). 
In dimcnsionlcss units, without loss of generality, we will 
assume that the blast wave propagates into the density 



Po(r) 



The shock radius is 



R = F 



(1) 



(2) 



The shock jump conditions give the density, velocity, and 
pressure at r = R — 0: p = 4po(R), v = (3/4)_R, p = 
(3/4:)po(R)R 2 . Then the density, velocity, and pressure 
of the downstream gas are of the form 



p = 4t-°0p(O 



(3) 



*Here and in the rest of the paper only the astronomically 
interesting case of adiabatic index 7 = 5/3 is considered. Gen- 
eralization to other values of 7 is straightforward. 



v = (3/4)^-^(0, 
p = (3/4)/3 2 t 2/3 - 2 - Q/5 p(0. 



(4) 



(5) 
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where p(l) = v (1) = p(\) = 1. 

When (3)-(5) are substituted into the gas dynamics 
equations (continuity, momentum, and pressure equa- 
tions), a system of ordinary differential equations for 
p(£)> an< i 1S obtained. With ' = d/d£, 



Then (6) gives 



(4f - 3v)p'/p = 3{v' + 2v/i) - 4a, 



(4£-3v)t/=p7p + 4(1 -!//?)«, 



(6) 



(7) 



(4f - 3u)j//p = 5(i/ + 2u/0 + 4(2 - 2//3 - a). (8) 

One can introduce new variables, reducing (6)- (8) to a 
system of two homogeneous equations [2] . We put 



(9) 



Here C is proportional to the sound speed, V(l) = 
C(l) = 1, and, after re-defining ' = £ d/d(, = d/dln£, 
we get 



V 



C' = 



F(V, C) 
H(V, C) ' 



G(V,<7) 



(4 - 3V)H(V, C) ' 



(10) 



(11) 



where 



F = - (15C 2 - (4 - 3F) 2 + 4(1 - /? _1 )(4 - 3V)) V 

+4 ((a- 2(1 -/T 1 ))*? 2 , (12) 

H(V,C) = 5C 2 - (4 -3F) 2 , (13) 

G(F,C) = C( J F+(6F-4/3- 1 )i7). (14) 

The system (10), (11) is integrated as follows: start at 
ln£ = with V = C = 1, and integrate back to ln£ = 
— oo. The integration can terminate at the singular point 
H = 0, unless the singularity is removed by having F = 
at this point, t 

Type I, [2]: For a < 3, energy conservation gives 
/? = 2/(5 — a). The blast wave is decelerating, /3 < 1. If 
also a > 2, the solution terminates at the singular point 
V = 4/3, C = at a finite value of ln£ = ln£ m i n , and 
an evacuated region forms at small £. Near the singular 
point, that is for small x = ln£ — ln£ m i n , equations (10), 
(11) give 



V = 4/3-(4/5)(a-l)x, Ctxi 1 / 2 . 



(15) 



^The singularity 4 — 3V = with H / does not occur for 
the values of a that we consider. 



p OC X 



8a - 18 
3a — 3 



(16) 



For a < 3, the density slope v < 1. The density singular- 
ity is intcgrable, and the pressure at the singular point 
p oc pC 2 is vanishing. 

Type II, [7]: For a > 3.26, (3 is determined numer- 
ically, by removing the singularity at H = 0. The self- 
similar solution extends all the way to ln£ = — oo. The 
blast wave is accelerating, (3 > 1. 

Type III, the new regime: For a > 3, the type I 
solution, (3 = 2/(5 — a) cannot be correct. This solution 
has infinite pressure at the termination point at finite £. 
Infinite pressure is not necessarily unphysical. We are 
considering just an asymptotic solution at infinite time. 
For a > 3 a finite fraction of the mass of the gas is not de- 
scribed by the self-similar solution, because it originates 
in the non-sclf-similar part of the unperturbed density 
(for a discussion sec [7]). In terms of pure self-similar 
solution, this gas forms an infinitely heavy piston, which 
might be able to support the infinite pressure at the end 
of the self-similar flow. However the problem is that the 
piston moves at a constant velocity, and therefore it can- 
not keep up with the accelerated blast wave. We must 
search among different values of (3. 

For a > 3, because of the constant-velocity piston 
formed by the non-self-similar gas, the solution cannot 
be decelerating, that is (3 > 1. It turns out (by numeri- 
cal integration of (10), (11)) that for a < 3.26 and (3 > 1 
the solution terminates at the singular point V = 4/3, 
C = at a finite value of ln£ = ln£ m ; n . Near the singu- 
lar point, that is for small x = ln£ — ln£ m j n , equations 
(10), (11) give for j3 > 1 



y = 4/3- (8/5)(2-/3- 1 )x, Cocx 1 / 2 . 



Then (6) gives 



p (X x u , v = 1 + 



5(a - 3) 
6(2-/3-!) 



(17) 



(18) 



For a > 3, the density slope v > 1. The pressure at the 
singular point p oc pC 2 is diverging. But this is impossi- 
ble for an accelerating solution, as discussed above. 

We are left with just one possibility, (3 = 1. In this 
case, (10), (11) give near the singular point 



V = 4/3 - 4(1 - a/5)x, 



Cocrr/ 1 , p 



Then (6) gives 



p oc x 



3(5 -a)' 



(19) 



(20) 



(21) 
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FIG. 1. The a = 2.85 density case. Thick lines show 
the self-similar solution with f3 = 2/(5 - a) = 0.930. Thin 
lines show the numerical solution. The horizontal axis gives 
?7-R s hock. Shown are 0.1p//9 s hock, l + «/fmax,2 + p/p max . Res- 
olution and run time: 10000 grid, 20 doublings of the grid. 
See main text for details. 



The pressure at the singular point p oc pC 2 is finite. 
Since the blast wave moves at a constant speed, this end 
pressure can be provided by the piston of non-self-similar 
gas. Since 2/i > 1, the mass is diverging, but this is 
OK, because in terms of an ideal self-similar solution the 
total mass of the exploding gas is infinite. Part of this 
mass forms an infinite-mass piston, another part forms 
an infinite-mass self-similar flow. 

The new solution exists only in a narrow interval of 
density slopes, 3 < a < 3.26. For a > 3.26, numerical 
integration of (10), (11) with (3=1 does not terminate 
at the point V = 4/3, C = 0. The solution crosses the 
line H = in a different point, with f ^ 0. Therefore, 
our solution is no longer valid. The right solution is the 
type II solution of [7]. 

Numerical simulation — There is no theorem guar- 
anteeing that a self-similar solution is asymptotically 
reached starting from a generic initial condition. One 
has to assume this, or check numerically. We numeri- 
cally simulated spherically symmetrical explosions in a 
power-law atmosphere. A spherically symmetric version 
of the code [9] was used. It was found that the second- 
order code becomes unstable near the termination part 
of the a w 3 blast waves. We therefore had to use the 
first-order version of the code. 

It was found that for a w 3 a very long run is needed 
to reach the theoretically predicted self-similar solution. 
This can be explained as follows. The energy conser- 
vation gives, approximately, M(R)R 2 — const, where 
M(R) is the enclosed mass 



FIG. 2. Same as fig. 1, for a = 3.15, (3 = 1. 



M(R) = J drr- a oc i? 3 ~ a - 1. (22) 

Here 1 is the length of the numerical grid, and all the 
energy of the explosion is put into the first grid point at 
the initial time. If a is close to 3, one has to wait for a 
long time, when the shock moves to large R, before one 
of the terms in (22) becomes much larger than the other. 
Since our new solution exists only in narrow interval of 
a, it is not expected to be very different from the Sedov- 
Taylor or Waxman-Shvarts solutions. To detect the small 
difference, we have to aim at a few percent accuracy for 
| a - 3| ~ 0.15. Putting i?l 3 ~ Q l = 30 gives the greed 
of ~ 10 10 points. This is more than we can simulate 
directly. We used much smaller grids, resolutions of 10 4 
and 10 5 , but then the grid was doubled for up to twenty 
times, when the blast wave was reaching the end of the 
grid. 

We first made sure that this numerical method does 
reproduce the hollow Sedov- Taylor solutions (2 < a < 3, 
(3 = 2/(5 — a)). It does, as figure 1 shows. We 
also checked that the numerical method reproduces the 
Waxman-Shvarts solution (it does for a = 3.40 with 20 
doublings of the grid) . 

We then checked if the new type III solution (3 < a < 
3.26, (3 = 1) is the true asymptotic of the explosion. It 
apparently is, as figure 2 shows. However, we have to 
say that our numerical proof is only tentative. It turns 
out that increased resolution makes the rear zone of the 
blast wave (numerically ?) unstable, fig. 3. The den- 
sity and velocity profiles do not change with increased 
resolution and still agree with the predicted self-similar 
profiles. The fluctuating pressure is still centered around 
the predicted self-similar solution. But the lack of nu- 
merical convergence near the end of the blast wave (or 
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FIG. 3. a = 3.15, analytical solution is not shown, 125000 
grid, 18 doublings of the grid. 



a possible indication of a true weak instability?) must 
somewhat reduce our confidence that our self-similar so- 
lution is the true long-time asymptotic of the explosion. 
Conclusion 

A new, type III, self-similar solution of the strong ex- 
plosion problem in the power-law atmosphere is given. 
This solution describes a constant speed ((3 = 1) blast 
wave, and applies in a narrow interval of density slopes, 
3 < a < 3.26. Numerical simulations confirm (tenta- 
tively) that the new solution is indeed a long-time asymp- 
totic of the strong explosion. 
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